#################
## Replication File for DeScioli and Kimbrough, "Alliance formation in a side-taking experiment", JEPS
#################
## import libraries
library(gtools)
library(foreign)

## create function to get 95% CI for later use
get.conf.int=function(x)t.test(x)$conf.int

## read in raw experimental data
x=read.csv("~/Dropbox (Personal)/Alliances/Experimental Paper/JEPS Final/AllSessions_JEPS.csv")
NUM_SESSIONS=length(unique(x$Session))

## create other variables for data set
id=rep(rep(1:8,each=21),NUM_SESSIONS)
session=rep(1:NUM_SESSIONS,each=8*21)
id=session*100+id
treatment=c(rep("Table",21*8*2),rep("NoTable",21*8*2),rep("Table",21*8),rep("NoTable",21*8),rep("Table",21*8),rep("NoTable",21*8),rep("Table",21*8*2),rep("NoTable",21*8*2),rep("ChatTable",21*8),rep("ChatNoTable",21*8*2),rep("ChatTable",21*8*2),rep("ChatNoTable",21*8),rep("ChatTable",21*8*3),rep("ChatNoTable",21*8*3),rep("ChatTable",21*8))
period=rep(1:21,8*NUM_SESSIONS)
full=ifelse(treatment=="Table",1,ifelse(treatment=="ChatTable",1,0))
chat=ifelse(treatment=="ChatNoTable",1,ifelse(treatment=="ChatTable",1,0))

## create data storage objects
sess_num=rep(1:NUM_SESSIONS,each=21)
supporterdiffarray=array(dim=c(168,8,NUM_SESSIONS))
lagalliancearray=array(dim=c(21,8,NUM_SESSIONS))
rankdiffarray=array(dim=c(21,7,NUM_SESSIONS))
lagbandwagonarray=array(dim=c(21,8,NUM_SESSIONS))
mutualarray=array(dim=c(21,7,NUM_SESSIONS))

## loop over the sessions to create measures of alliance structures
for(j in 1:NUM_SESSIONS){
	fightMat=diag(0,8,8)
	fightMat2=NULL
	lagallianceMat=matrix(nrow=21,ncol=8)
	rankdiffMat=matrix(nrow=21,ncol=7)
	lagbandwagonMat=matrix(nrow=21,ncol=8)
	mutualMat=matrix(nrow=21,ncol=7)

	for(i in 1:21){
		y=x[intersect(which(sess_num==j),which(x$Period==i)),]
		
		## transform rank vector into a matrix
		rankMat1=matrix(as.numeric(y[8:63]), nrow=8,ncol=7, byrow=T)
		rankMat=diag(0,8,8)
		for(f in 1:8){
			if(f==1)rankMat[f,]=c(0,rankMat1[1,])
			if(f>1&&f<8)rankMat[f,]=c(rankMat1[f,1:(f-1)],0,rankMat1[f,f:7])
			if(f==8)rankMat[f,]=c(rankMat1[8,],0)
		}
		
		## construct measure of rank differences (asymmetry between i and js ranks of one another)
		rankDiffs=matrix(nrow=8,ncol=7)
		for(q in 1:8){
			for(r in 1:7){
				myranker=which(rankMat[q,]==r)
				if(length(myranker)==1){
					yourrank=rankMat[myranker,q]
					rankDiffs[q,r]=abs(yourrank-r)
					}
			}
		}
		rankdiffMat[i,]=colMeans(rankDiffs,na.rm=T)
		
		## create indicator matrix with 1s when i and j's ranks are mutual
		mutualMat[i,]=colSums(matrix(as.numeric(rankDiffs==0),nrow=8),na.rm=T)
		
		if(i > 1){
			y=x[intersect(which(sess_num==j),which(x$Period==(i-1))),]
			rankMat2=matrix(as.numeric(y[8:63]),nrow=8,ncol=7,byrow=T)
			rankMat_prev=diag(0,8,8)
			for(f in 1:8){
				if(f==1)rankMat_prev[f,]=c(0,rankMat2[1,])
				if(f>1&&f<8)rankMat_prev[f,]=c(rankMat2[f,1:(f-1)],0,rankMat2[f,f:7])
				if(f==8)rankMat_prev[f,]=c(rankMat2[8,],0)
			}
			cors=rep(0,8)
			
			## compute alliance building measure
			for(g in 1:8){
				cors[g]=cor(rankMat[g,][-g],rankMat_prev[,g][-g])
				if(is.na(cors[g]))cors[g]=0
				
							}
			lagallianceMat[i,]=cors
			
		}
		## determine players' strengths, i.e the share of players that would support i in each possible fight against each j
		for(k in 1:7){
			for(l in (k+1):8){
				favork=1
				favorl=1
				fighters=c(k,l)
				all=setdiff(1:8,fighters)
				for(m in all){
					if(rankMat[m,l]<rankMat[m,k]){
						favorl = favorl + 1
						}
					else{
						favork = favork + 1
						}
					}
				fightMat[k,l]=favork-favorl
				}
			}
		fightMat2=rbind(fightMat2,fightMat)
		
		## compute measure of bandwagon strategy (correlation of i's ranks at t with others' ranks of j at t-1)
		if(i>1){
			for(q in 1:8){
				cors[q]=cor(rankMat[q,][-q],(colSums(rankMat_prev)-rankMat_prev[q,])[-q])
				if(is.na(cors[q]))cors[q]=0
				}
			lagbandwagonMat[i,]=cors	
		}
		
		
	}
	## store data in arrays
	supporterdiffarray[,,j]=fightMat2
	lagalliancearray[,,j]=lagallianceMat
	lagbandwagonarray[,,j]=lagbandwagonMat
	rankdiffarray[,,j]=rankdiffMat
	mutualarray[,,j]=mutualMat
}

ind_data=as.data.frame(cbind(id, session,period= c(period),  alliance=as.vector(lagalliancearray),egalitarian_bandwagon=as.vector(lagbandwagonarray),treatment=as.factor(treatment), chat, full))

########################
## Results
########################

########################
## Alliance, bandwagon and egalitarian strategies in periods 11-20
########################

### NOTE:
## Chat/Partial Info treatment = 1, Chat/Full Info treatment = 2, No Chat/Partial Info = 3, No Chat/Full Info = 4

### perform analyses on individual-level means
## get subject level means periods 11-20 (coded as 12-21 in the data)
individual_means=aggregate(cbind(alliance,egalitarian_bandwagon)~id+treatment+chat+full,FUN="mean", data=ind_data, subset=c(period>11))

## t-tests for alliance building
for(i in c(4,2,3,1)){print(t.test(individual_means$alliance[which(individual_means$treatment==i)]))}

## Anova for alliance building
reg1=lm(alliance~chat*full,data=individual_means)
anova(reg1)

## t-tests for egalitarian/bandwagon
for(i in c(4,3,2,1)){print(t.test(individual_means $egalitarian_bandwagon[which(individual_means$treatment==i)]))}

## Anova for egalitarian/bandwagon
reg2=lm(egalitarian_bandwagon~chat*full,data= individual_means)
anova(reg2)

## Create data for Figure 1
treatment_means=aggregate(cbind(alliance,egalitarian_bandwagon)~chat+full,FUN="mean", data=ind_data, subset=c(period>11))
treatment_ci=aggregate(cbind(alliance,egalitarian_bandwagon)~chat+full,FUN="get.conf.int", data=ind_data, subset=c(period>11))

treatment_means
treatment_ci

#######################
## Side Taking Dynamics
#######################
##### Create data for Figure 2
treatment_means_by_period=aggregate(cbind(alliance,egalitarian_bandwagon)~period+chat+full,FUN="mean", data=ind_data)

treatment_means_by_period

###########
## TABLE 1 and TABLE A1
###########
## many regression analyses are more automatic in Stata, so we conduct the analysis underlying Table 1 using Stata
## create type classifications for analysis in Appendix A, table A1
type_classification=rep(NA, nrow(ind_data))
for(i in 1:nrow(ind_data)){type_classification[i]=which.max(c(ind_data$alliance[i], ind_data$egalitarian_bandwagon[i], abs(ind_data$egalitarian_bandwagon[i]), 0.2))}

## store data
ind_data=cbind(ind_data,type=type_classification)

## save data for Stata regressions
write.dta(ind_data,"~/Dropbox (Personal)/Alliances/Experimental Paper/JEPS Final/Individual_Data_JEPS.dta")

## see other analysis file for Stata regressions

########################
## Best Friends and Enemies
########################
treatment=c(rep("Table",21*2),rep("NoTable",21*2),rep("Table",21),rep("NoTable",21),rep("Table",21),rep("NoTable",21),rep("Table",21*2),rep("NoTable",21*2),rep("ChatTable",21),rep("ChatNoTable",21*2),rep("ChatTable",21*2),rep("ChatNoTable",21),rep("ChatTable",21*3),rep("ChatNoTable",21*3),rep("ChatTable",21))
## create a data frame with the percentage of mutual ranks
mat=NULL
for(i in 1:NUM_SESSIONS){
	mat=rbind(mat,mutualarray[,,i])
}
mutual=as.data.frame(cbind(period=rep(1:21, NUM_SESSIONS),as.factor(treatment),mat))
colnames(mutual)=c("period","treatment","first","second","third","fourth","fifth","sixth","seventh")

## Chat/Partial Info treatment = 1, Chat/Full Info treatment = 2, No Chat/Partial Info = 3, No Chat/Full Info = 4

## compute data underlying figure 3
percent_mutual=aggregate(cbind(first,second, third,fourth,fifth,sixth,seventh)~treatment,data=mutual,FUN="mean", subset=c(period>11))/8

percent_mutual

## compute Chi-Sq tests
mutual_counts=aggregate(cbind(first,second, third,fourth,fifth,sixth,seventh)~treatment,data=mutual,FUN="sum", subset=c(period>11))
## test of equal proportions with each rank 
for(i in 1:4)print(chisq.test(rbind(mutual_counts[i,2:8],560-mutual_counts[i,2:8])))

## Chat partial vs chat full
prop.test(x=c(sum(mutual_counts[1,2:8]), sum(mutual_counts[2,2:8])), n=rep(560*7,2))

## no chat partial vs no chat full
prop.test(x=c(sum(mutual_counts[3,2:8]), sum(mutual_counts[4,2:8])), n=rep(560*7,2))

## comparing rank 7 in partial vs full
mutual_counts

## chat, populated with data from the mutual_counts table
prop.test(x=c(73,162), n=c(560,560))

## no chat, populated with data from the mutual_counts table
prop.test(x=c(68,144), n=c(560,560))


########################
## Chat Messages and their Consequences
########################
## Power Differences

## create measure of average power difference
average_power_difference=rep(NA,NUM_SESSIONS*21)
count=0
for(j in 1:NUM_SESSIONS){
	for(i in 1:21){
		count=count+1
		starts=seq(1,168,8)
		mat=supporterdiffarray[starts[i]:(starts[i]+7),,j]
		## average power difference = average asymmetry of alliance size over all 28 pairwise fights
		average_power_difference[count]=sum(abs(mat))/((8*(7))/2)
	}
}

full=ifelse(treatment=="Table",1,ifelse(treatment=="ChatTable",1,0))
chat=ifelse(treatment=="ChatNoTable",1,ifelse(treatment=="ChatTable",1,0))

## store as data frame
powerdiffs=as.data.frame(cbind(period=rep(1:21, NUM_SESSIONS),as.factor(treatment),average_power_difference, full, chat, session=rep(1:NUM_SESSIONS,each=21)))

## Produce data for Figure 4
pdiff_mean=aggregate(average_power_difference~full+chat+session,data=powerdiffs,FUN="mean", subset=c(period>12))
pdiff_ci=aggregate(average_power_difference~full+chat+session,data=powerdiffs,FUN="get.conf.int", subset=c(period>12))

pdiff_mean
pdiff_ci

